The locomotion of living animals, such as flying birds and swimming fishes, has persistently been the focus of a multidisciplinary research effort, from ethologists and biologists to physicists and applied mathematicians. The previous studies paid much attention to explore the optimal propulsive force and efficiency generated by a slender body with prescribed oscillation,^[1–3] and the experimental and numerical results on the flapping foil were in good agreement with those from the biological research on living aquatic and flying animals.^[4,5] Thus, it is reasonable to interpret that the hydrodynamic or aerodynamic performance may play a vital role in the natural selection.

During the flapping motion, the foil would experience the resistive force (drag) or the propulsive force (thrust) which is mainly decided by the Strouhal number St = 2fA/U, where f is the flapping frequency of the foil, 2A denotes the width of the wake behind the foil, and U represents the velocity of swimming speed. In the Newtonian fluid, the drag is always related to the Bénard–von Kármán wake produced by the flapping foil while the reverse Bénard–von Kármán (rBvK) wake leads to the thrust on the foil.^[6] For an ideal flow, the transition between the drag to the thrust is synchronous with the BvK–rBvK transition when the St exceeds an intermediate value.^[7] However, Godoy-Diana et al.^[8] experimentally uncovered that the asynchrony exists in the force and wake transitions of a pitching foil in the viscous flow. Their results showed that the transition from a BvK wake to a rBvK wake precedes the actual drag–thrust transition, and the rBvK wake corresponds to the resistive force at lower St, which is contradictory to the classical inviscid theory. Moreover, the authors also found the symmetry breaking of the rBvK wake at higher St which results in an asymmetric pattern. Those conclusions have been validated by other researchers.^[9–12] Deng and Caulfield,^[13] and Deng et al.^[14,15] expanded the process of wake transition based on the previous studies, they observed that a two-dimensional (2D) deviated wake behind the 2D foil would show the features of the three-dimensional (3D) wake when the St keeps increasing.

Although the asynchrony of the BvK–rBvK transition and drag–thrust transition shows an exciting phenomenon corresponding to the difference between the in-viscid and viscous fluid, researchers may focus more on the stability of the wake produced by the flapping foil, in the other word, the process of the transition from a symmetric rBvK wake to the asymmetrical wake. To the best of our knowledge, Bratt^[16] firstly found the asymmetry pattern in the wake of a pitching foil. Depending on the experimental and numerical studies,^[17] Jones et al. also showed the deviation of the symmetry wake occured at a plunging foil. The determining factor on the deviated direction of the asymmetrical wake was also investigated by Lai and Platzer,^[18] von Ellenrieder and Pothos,^[19] and Zheng and Wei.^[20] The previous works show that the deviation direction depends strictly on the initial conditions of the plunging motion, the wake would deviate upwards when the foil plunged to the downwards at the beginning of the oscillation, and vice versa. Godoy-Diana et al.^[21] proposed a theoretical model to predict the wake deviation based on the relationship between the real phase velocity of the vortex and the self-advection velocity of two consecutive counter-rotating vortices, which allowed Jallas and Marquet^[22] to analyze the linear and nonlinear perturbations on the wake symmetry-breaking behind a pitching foil. Moreover, the studies on a flexible foil showed that the low-flexibility foil would maintain the stability of the wake structures^[23] while the high-flexibility one instead enhances the deviation of the wake.^[24]

Considering the literatures mentioned above, we find that the all studies focused on the wake structures of a pitching foil were conducted on the symmetry foil. To our best acknowledge, however, there are limited studies paid attention to the effect of foil’s shape on its wake structures. Since the asymmetric geometry of the foil is widely observed at the real animals, we use a 2D model to reveal the inherent relationship among the shape’s asymmetry, force generation, and wake structures of a pitching foil. When the aspect-ratio (AR) of a foil is larger enough, the 2D model can uncover the key dynamic elements of the creation and organization of the vortex generated by a 3D flapping foil. Particularly, the previous work shows that the rBvK wake produced by a 2D foil is similar to the vorticity in the center plane along the span-wise direction of a 3D foil at AR ≥ 4.^[25,26] Therefore, the present 2D model can reveal the crucial information of a 3D pitching foil with larger AR since the corresponding 3D wake is quasi2D wake.

In our work, we give the problem description and numerical methods in Section 2, where the schematic of the flow configuration, kinematics of the oscillating foil, relevant parameters, and numerical approach are described. Then, we discuss the numerical results in Section 3. The last section is the closure of the paper in which we summarize the essential observations made from our study.

In this section, we merely introduce the present problem and the corresponding computational approach, the work of Godoy-Diana et al.^[8] is considered as the baseline of our study.

2.1. Problem description

The uniform flow around a 2D foil immersed in a fluid of kinematic viscosity ν with an oncoming velocity U is investigated. The foil has a semicircle leading edge and a cuspidal trailing edge, as shown in Fig. 1.

Fig. 1.

	Figure Option ViewDownloadNew Window
	Fig. 1. Flow configuration and the asymmetrical foil: (a) sketch of the flow configuration, (b) design of the asymmetrical foil.

Similarity with Godoy-Diana et al.’s experiments,^[8] the thickness of the pitching foil calculated as D/c = 5/23 and the chord-based Reynolds number Re_D = DU/ν = 255 are used for the entire study. The center of the prescribed oscillated on the foil is located at D/2 from the leading edge, and we argue that the flapping motion follows the sinusoidal law 1 where θ₀ is the maximum rotational angle, and t denotes the instantaneous time. We also define two dimensionless parameters A_D and Sr to describe the kinematics of the pitching foil. The dimensionless peak-to-peak pitching amplitude A_D is 2 where A can be calculated as A = sin θ₀(c − D/2).

The chord-based Strouhal number Sr, also known as Keulegan–Carpenter number, is defined as 3 thus it is easy to understand that St is the product of Sr and A_D.

For the asymmetrical foil, we use the symbol K to describe the asymmetry of the geometry, which is defined as the ratio of the ordinate values of close point B₀ to the half value of the thickness. Therefore, as shown in Fig. 1(b), it is easy to find that the K is ranged from 1 (upper one) to −1 (bottom one) and K = 0 represents the symmetric foil (middle one).

2.2. Numerical method

Based on the previous studies,^[27] the nonlinear Navier–Stokes equations can be translated to a non-inertial frame of reference (e_X, e_Y), as illustrated in Fig. 1(a). Around the Z-axis satisfied the righthand rule on the XoY plane, thus the rotational speed is ω = dθ/dt but would not affects the frame of reference (e_x, e_y). The 2D incompressible Navier–Stokes (NS) equations are described as 4 5 where u = (u, v)^T denotes the vector field, p is the pressure, ω (t)e_z × u represents the rotational acceleration, and the convective velocity w is defined at spatial coordinates X = (X, Y) as 6 where the first term represents the translational velocity and the second term denotes the angular velocity of the oscillated foil. Therefore the boundary condition on the foil surface is 7

In the present work, the computational domain is considered as a rectangular area. The inlet velocity boundary is at 8c from the leading edge of the foil, as u = (U, 0)^T and ∂p / ∂x = 0. The is used as the nostress outflow boundary condition, and the downstream pressure outlet is at 20c from the trailing edge of the foil. For the upper and lower flow boundaries located at 8c from the certain line of the foil, the slip wall condition is used.

The finite element method is used to achieve the spatial discretization of the N–S equations (4) and (5), and the temporal discretization is executed with the second-order accuracy backward-implicit scheme. The SIMPLE algorithm is used to achieve the pressure–velocity coupling of the continuity equation (4). Moreover, Gauss–Seidel linear equation solver is employed to solve the mentioned discretized equations. This numerical method, implemented in the available commercial software FLUENT 14.0, is used to simulate the present problems, and one can find the validated process of this method in our published papers.^[28]

The prescribe motion equation (1) is handled via the source term equaled to the angular acceleration of the pitching foil, and the source term is added in the flow governing equations using a user defined function (UDF). The DEFINE_CG_MOTION function is used to smooth and regenerate the cells around the foil at each updated time, and the time step is set to 0.0001c/U in per period to ensure the grid quality at each computational step. To save the computational time, the adjacent area of the foil is surrounded by the unstructured mesh while the rest of the computational domain is set as the structured mesh. Moreover, the adjoining grid of the body is distributed with higher resolutions to precisely capture the hydrodynamic force and the wake structures of the pitching foil. The first cell is put away from the foil surface at a reasonable distance of 0.0001c with y⁺ = O(1).

To validate the present numerical predictions, we compare the thrust generation and wake structures produced by a symmetric foil with the previous studies performed by Godoy-Diana et al.^[8] As shown in Fig. 2(a), the result shows that there are no significant differences between present and the published paper. The grid independence tests are performed with three mesh sizes: a coarse mesh of 1.0 × 10⁵, a medium mesh of 4.0 × 10⁵, and a fine mesh of 4.8 × 10⁵, corresponding to 201, 421, and 501 points on the surface of the foil, respectively. The results of the tests indicate that the medium mesh is most suitable to our study because it obtains higher accuracy with less number of grids.

Fig. 2.

	Figure Option ViewDownloadNew Window
	Fig. 2. (color online) (a) Validation of present numerical approach: (a1) the present numerical result, (a2) reproduced from Godoy-Diana et al.,[8] copyright 2008, with permission from APS. The symbols used in (a) include □ : BvK wake; +: rBvK wake; ■ : aligned vortices; △: asymmetric wake. (b) The grid independence tests.

To investigate the effect of the asymmetrical geometry on the wake structures of a pitching foil, the nonlinear simulations are performed at the Sr ranging from 0.1 to 0.5, 0.25 ≤ A_D ≤ 2.5, and 1 ≤ K ≤ 1 with the interval 0.1. Here, we only show the numerical results for K = 1 and K = 1, but it is noted that all simulations reveal similar changing rules of the wake structures among different values of K.

The evolution of wake structures, evidently shown in the (Sr, A_D) plane (Fig. 3), illustrates a similar process for K = 1 and K = 1. The P+S wake (symbols + in Fig. 3) denoted one pair and one single vortex shed by the foil in one period^[29] is initially observed at the lower left side in the (Sr, A_D) plane. It is easily found that the P+S wake always exists at lower values of A_D and Sr, such as A_D = 0.25 or Sr = 0.1, which agrees well with previous studies.^[30] After slight increase of A_D and Sr, the P+S wake evolutes to the orthodox BvK wake (□ symbols) in which the vortices are shed from alternating sides of the object, and the clockwise vortices (blue ones in the vortices map) move above the horizontal line and the counterclockwise ones (red ones in the vortices map) move below. Then, the positive (clockwise) and negative (counterclockwise) vortices are aligned with the symmetry line of the wake, as a 2S wake (■ symbols). The 2S wake is also considered as the transition boundary (red dash-dot line in Fig. 3) between the Bvk wake and the rBvK wake. By the sequential increase of A_D and Sr, the vortices shift to an opposite rotational direction and transit into the rBvK (○ symbols). Therefore, the mean flow behind the foil has a typical jet profile corresponding to the thrust.^[31] When further increasing the Sr and A_D, the stability of the symmetric rBvK wake is broken, and the deviation of the wake occurs as an asymmetrical wake (▽ symbols).

Fig. 3.

Figure Option ViewDownloadNew Window

Fig. 3. (color online) Two examples of the numerical points in the (Sr, AD) plane: (a) the wake phase diagram for K = 1, (b) the wake diagram for K = 1. +: P+S wake; □: BvK wake; ◯: rBvK wake; ▽: asymmetric wake; ■: aligned vortices (2S wake). Red dashdot line: BvKr–BvK transition boundary; black solid line: transition between the rBvK wake and the asymmetrical wake for an asymmetrical foil; black dash line: transition between the rBvK wake and the asymmetrical wake for a symmetrical foil (K = 0).

It is illustrated that the symmetric-breaking (black solid line in Fig. 3) of the rBvK pattern for an asymmetrical foil is not identical with the transition boundary (black dash line in Fig. 3) between the rBvk wake and the asymmetrical wake for a symmetrical foil. In general, one can observe that the transition from the rBvK wake to the asymmetrical wake for an asymmetrical foil precedes the transition for a symmetric foil when K < 0, and vice versa. As the sign of K describes the asymmetrical direction of the foil, we argue that the wake deviation is handled by the direction of the initial motion and the asymmetrical direction of the foil. When those directions are the same, the wake structures behind a pitching foil are easily to be instable and result in the asymmetrical wake. However, the asymmetrical geometry of the foil would refrain from the stability-breaking of the rBvK wake when the mentioned directions are opposite. From the perspective of the vortex generation, the asymmetric geometry of the foil leads to the different initial velocities of the vortex shedding from the foil’s up/down side. Therefore, the induced velocities between two vortices produced in one period are not equal, and this situation would advance/delay the deviation of the symmetric rBvK wake. When the net induced velocity has the same direction with the deviated direction of the wake, the symmetry-breaking would be advanced, and vice versa. This finding may be used to understand that the living beings change their propulsive organism to an asymmetrical geometry to actively control the flow field.^[32,33]

Moreover, we theoretically consider the effect of asymmetrical geometry on the wake’s instability. For the arbitrary asymmetrical foil, a predicted work is performed to access the deviation of the rBvK wake in the incompressible Newtonian fluid with a uniform laminar flow. As shown in Fig. 4, an asymmetrical foil with special K leads a θ₁ between the upper boundary of the foil and the horizontal axis, as well as a θ₂ between the bottom boundary of the foil and the horizontal axis. Thus, we have tan θ₁/tan θ₂ = K(1 − K). The difference between θ₁ and θ₂ can be used to estimate the positive or negative value of K, such as θ₁ − θ₂ > 0 describes K < 0, θ₁ − θ₂ < 0 indicates K > 0, and θ₁ − θ₂ = 0 denotes a symmetrical foil.

Fig. 4.

	Figure Option ViewDownloadNew Window
	Fig. 4. (color online) A sketch of symmetric-breaking for rBvK.

We now consider the critical state of the symmetry-breaking, where a rBvK wake produced by the pitching foil is propagating to the downstream. The clockwise-rotation vortex is seemed like the positive vortex with circulation Γ (red ones), while the circulation of the anticlockwise one is Γ (blue ones). The b is the vertical gap between the consecutive vortices, while a is the horizontal gap and ξ is the gap between the vortex cores. Considering those two point vortexes as a dipolar structure, we have the following equation:^[34] 8 We define the initial speed of the vortices shed from the foil as U_ini, and recall the work of Godoy-Diana et al.^[21] on the wake instability of the symmetrical foil. Thus the criterion of the instability of the wake produced by the asymmetrical foil is 9 Here the first two terms are related to the symmetric-breaking for the symmetrical foil, and the last term corresponds to the asymmetrical geometry effect. With this equation, we show that the symmetric-breaking for the asymmetrical foil occurs at . Moreover, equation (9) indicates how the deviation of the wake is handled by the kinematics of the foil (U_dip, α, and U_ini), the flow conditions (U), and the geometry (θ₁ and θ₂) of the foil. When θ₁ > θ₂, for example, it is ready to obtain , so that the deviation of the wake produced by asymmetrical foil K < 0 proceeds the symmetric-breaking of the rBvK wake generated by the symmetrical foil. On the other hand, for θ₁ < θ₂, we obtain the opposite results that the symmetric-breaking of the rBvK wake is lagged.

With the discussion above, we can predict the symmetry-breaking of a rBvK wake for the asymmetric foil through solving the present theoretical formulation. Particularly, our work is based on the conclusions that the wake generated by a 3D foil is quasi2D at larger aspect-ratio. Therefore, the present model can be well-used when AR ≥ 4.

We numerically solve the nonlinear Navier–Stokes equations to understand the effect of asymmetrical geometry on the wake structures of a pitching foil. Some typical wake structures are plotted in the (Sr, A_D) plane to show the relationship among the kinematics of the foil, the foil’s asymmetry, and the fluid dynamics around the foil. The numerical results uncover that the asymmetric geometry of foil’s shape plays an important role affecting the symmetry-breaking of the rBvK wake generated by the foil. The rBvK wake is easier to become instable when the asymmetrical direction of the foil and the pitching direction of the initial motion are the same, while the deviation of the symmetric wake is delayed when the two directions are opposite. Moreover, we formulate a theoretical model to predict the symmetry-breaking and provide some helpful findings for designing the bio-inspired robotics.